Second order conformal multi-symplectic method for the damped Korteweg–de Vries equation
Guo Feng
Fujian Province University Key Laboratory of Computational Science, School of Mathematical Sciences, Huaqiao University, Quanzhou 362021, China

 

† Corresponding author. E-mail: hydhgf@163.com

Project supported by the Program for Innovative Research Team in Science and Technology in Fujian Province University, China, the Quanzhou High Level Talents Support Plan, China (Grant No. 2017ZT012), and the Promotion Program for Young and Middle-Aged Teacher in Science and Technology Research of Huaqiao University, China (Grant No. ZQN-YX502).

Abstract
Abstract

A conformal multi-symplectic method has been proposed for the damped Korteweg–de Vries (DKdV) equation, which is based on the conformal multi-symplectic structure. By using the Strang-splitting method and the Preissmann box scheme, we obtain a conformal multi-symplectic scheme for multi-symplectic partial differential equations (PDEs) with added dissipation. Applying it to the DKdV equation, we construct a conformal multi-symplectic algorithm for it, which is of second order accuracy in time. Numerical experiments demonstrate that the proposed method not only preserves the dissipation rate of mass exactly with periodic boundary conditions, but also has excellent long-time numerical behavior.

1. Introduction

The Korteweg–de Vries (KdV) equation

was initially derived as a model of unidirectional propagation of water waves with small amplitude in a channel. Since the work of Korteweg and de Vries, it has been used in a large variety of physical situations, such as an acoustic wave in an anharmonic crystal, magnetohydrodynamic waves in warm plasmas, and ion acoustic waves. In many real situations, however, one cannot neglect the dissipation of energy. To account for this effect, here, we consider the damped Korteweg–de Vries (DKdV) equation with initial boundary value conditions as follows:
where constant β is a dispersion coefficient, is a damping coefficient, and is a given smooth function.

It is easy to verify that equation (2) has the following energy balance equation:[1,2]

where the momentum (energy) . Using the above relations and the exact solution of Eq. (1), references [14] derived several approximate solutions for Eq. (2). Song et al.[5] provided some numerical and theoretical analyses for it. Also, the conservation of mass for Eq. (2) can be given by

Various numerical methods have been proposed to solve the KdV equation, including the finite difference methods,[6,7] the finite element methods,[8] the spectral and pseudo-spectral methods,[9] etc. In particular, the KdV equation is also studied with some structure-preserving algorithms, such as symplectic and multi-symplectic methods,[1015] which can not only preserve some physical properties of the original problem under appropriate discretizations,[1618] but also have long time numerical behavior. Nevertheless, conventional structure-preserving methods mainly focus on the systems without dissipation, and in the presence of damping or diffusion, those structures would break down. Recently, McLachlan and Perlmutter,[19] as well as McLachlan and Quispel[20] introduced a concept of conformal symplectic methods for linearly damped Hamiltonian ordinary differential equations (ODEs). Then, Moore et al.[21,22] generalized it to multi-symplectic partial differential equations (PDEs) with linear dissipation and proposed some conformal multi-symplectic methods accordingly. Such methods can well preserve the conformal properties, such as conformal symplecticity or multi-symplecticity, conformal energy or momentum, and other quantities associated with the dissipative systems. Compared with standard structure-preserving methods, conformal structure-preserving methods preserve the dissipation rate exactly. In this paper, by employing the Strang-splitting technique and Preissmann box scheme to the DKdV Eq. (2), we obtain a second order conformal multi-symplectic scheme for it, which exhibits excellent performance in preserving the dissipative property of the system.

The rest of the paper is organized as follows. In Section 2, we present the damped multi-symplectic formulation and the corresponding conformal conservation laws for the DKdV equation. In Section 3, some operators and their properties are given. And then, we construct a conformal multi-symplectic method for the DKdV equation and prove some conservative properties for the proposed method. In Section 4, numerical experiments are carried out to show the effectiveness and advantages of this method. Finally, we make some conclusions in Section 5.

2. Conformal multi-symplectic structure and conformal conservation laws for the DKdV equation
2.1. Damped multi-symplectic Hamiltonian system and conformal conservation laws

Many weakly dissipative Hamiltonian PDEs can be written as a damped multi-symplectic Hamiltonian system[21]

where and are constant skew-symmetric matrices, , is a smooth function of . The corresponding conformal multi-symplectic conservation law[22] is
where Methods that preserve Eq. (6) are called conformal multi–symplectic integrators. If a = 0, the corresponding conformal energy conservation law is
If b = 0, the corresponding conformal momentum conservation law is
Note that integrating Eq. (8) in space with periodic boundary conditions yields the total conformal momentum
where . A numerical method is said to preserve the total conformal momentum if it satisfies the discrete versions of Eq. (9).

2.2. Conformal multi-symplectic structure for the DKdV equation

By introducing auxiliary variables , we can reform the DKdV Eq. (2) in the conformal multi-symplectic PDEs as

where , and
The corresponding conformal multi-symplectic conservation law is
and the conformal momentum conservation law is
Integrating the conservation law in Eqs. (11), (12), and the first equation of Eq. (10) with respect to x, respectively, under the proper boundary conditions such as periodic boundary conditions, gives the following total conformal conservation laws:

3. Conformal multi-symplectic scheme for the DKdV equation
3.1. Operator definitions and properties

In order to conveniently derive the algorithms, we first introduce some notations: the uniform grids and , where and are spatial length and temporal step span, respectively. The approximation of the value of function at the node is denoted by . We also define some finite difference operators and averaging operators

They satisfy the following commutative laws:[23]

3.2. Splitting method

In order to construct a conformal multi-symplectic scheme, we split the damped multi-symplectic Hamiltonian system of Eq. (10) into a conservative Hamiltonian system

and a dissipative system
For the former system, we solve it with a standard geometric integrator, such as the Preissmann box scheme, which is stated by . While the dissipative part is solved exactly by the analytical flow . Thus, the damped Hamiltonian system of Eq. (10) can be solved by composing these two flow maps, e.g., the Strang method[24] . McLachlan and Quispel[20,25] have shown that splitting methods may be used to numerically preserve conformal symplecticity of conformal Hamiltonian systems. This idea was extended to conformal multi-symplectic PDEs in Ref. [26], using symplectic Euler discretizations; other similar splitting techniques have been employed in the context of multi-symplectic PDEs by Kong et al.[27,28] In our work, we use the Strang-splitting technique to solve the DKdV system of Eq. (10). Following the above illustration, we denote the solutions of every step in the Strang splitting as follows:
Then a conformal multi-symplectic method based on the Preissmann box scheme is obtained.

3.3. Conformal multi-symplectic method for the DKdV equation

In this subsection, we present a conformal multi-symplectic (CMS) method for Eq. (10)

Theorem 1 The discretization equation (18) admits the following conformal multi-symplectic conservation law:

where .

Proof Taking the wedge product of the variational formulation of Eq. (18) with yields

The first term on the left side becomes
Similarly, we have . So equation (20) can be reduced to . Multiplying throughout it by completes the proof.

Applying the conformal multi-symplectic method of Eq. (18) for the damped multi-symplectic structure of Eq. (10), we obtain the following scheme:

Eliminating the auxiliary variables gives an equivalent form as

According to Theorem 1, we arrive at Theorem 2.

Theorem 2 The discretization equation (21) or (22) is a conformal multi-symplectic scheme and satisfies the discrete conformal multi-symplectic conservation law

where .

Under the periodic boundary conditions, summing the first equation in scheme of Eq. (21) over all spatial index j yields the discrete total conformal mass

namely, the conformal multi-symplectic scheme of Eq. (22) preserves the dissipation rate of mass exactly. Here, we use to measure the error of the dissipation rate of mass. We also define the discrete total conformal momentum as

4. Numerical experiments

In this section, we conduct some numerical experiments to exhibit the performance of the proposed method in Eq. (22). All solutions are computed with periodic boundary conditions in [−50,50] and the time interval is taken as [0,60] for all experiments. The contents include: (i) to display the accuracy of solution; (ii) to simulate the propagation of one soliton and double solitons collision over a long time, respectively; (iii) to show the variation of global conformal mass and momentum and monitor the preservation of the dissipation rate of mass against time; (iv) to compare the error of dissipation rate of mass with the Preissmann box (PBX) scheme, i.e., .

The DKdV equation (2) possesses the following approximate soliton solution:

where k is initial amplitude. It represents the motion of the single solitary wave with amplitude , moving with velocity and spatial width of the soliton . We take
where c is the integration constant, which is related to the initial location of the solitary wave.

The temporal accuracy of the proposed method of Eq. (22) is tested with[29]

where and denotes the numerical result at T with time step size τ. is L2 or L norm. Table 1 exhibits the computational results at T = 10 corresponding to a=0.01, 0.03 with the fixed space step length h=0.01 and different time step sizes. The parameters are fixed by β=0.1, k=0.3, c = 1. As is expected, the convergence order is consistent with our academic analysis of .

Table 1.

Temporal convergence order with h=0.01 by CMS in Eq. (22).

.

Similarly, we test the spatial accuracy with . Table 2 presents the computational results at T = 10 with various spatial steps and fixed temporal size τ=0.01. Other parameters are the same as before. It suggests that the spatial approximation order is close to 2.

Table 2.

Spatial convergence order with τ=0.01 by CMS in Eq. (22).

.

Example 1 Single soliton. Many physical systems, which are nonlinear, dissipative, and dispersive, can be modeled successfully using the damped KdV equation. Here, we consider a damping quantum fluid model. We would like to investigate the damping and dispersion effects on the formation of the solution. In order to observe the damping effect, we consider the initial condition

where h=0.1, τ=0.01, and c = 1. Figure 1 and 2 display the evolution of single solitons with damping coefficient a=0.01 and a = 0.02, respectively. The comparison of CMS and PBX in preserving the dissipation rate of mass is also presented. As we can see, the propagation of single solitons presents a good preservation of the phase space structure. From Figs. 1(b), 1(d), 1(e), 2(b), 2(d), and 2(e), it is obvious that the damped soliton moves with a gradually diminishing velocity and the soliton amplitude decreases over time. This implies that there exists kinetic energy loss because of the damping effect. In addition, the bigger the damping coefficient is, the faster the amplitude of the soliton decays, and the slower the soliton propagates, which means that the kinetic energy loss increases when the damping effect increases. Furthermore, the total conformal momentum and mass both decrease exponentially with time. Figure 1(c), 1(f), 2(c), and 2(f) indicate that the CMS method can preserve the dissipation rate of mass precisely, while the PBX method cannot reach.

Fig. 1. Numerical simulation of Eq. (22) with a=0.01 and the comparison of methods CMS and PBX with a=0.1. (a) Propagation of single solitary wave u; (b) numerical solutions at t = 0, t = 30, and t = 60; (c) error of the dissipation rate of mass r; (d) the evolution of total conformal momentum (e) the evolution of total conformal mass (f) comparison of methods CMS and PBX in the error of dissipation rate of mass with a=0.1.
Fig. 2. Numerical simulation of Eq. (22) with a=0.02 and the comparison of methods CMS and PBX with a=0.2. (a) Propagation of single solitary wave u; (b) numerical solutions at t = 0, t = 30, and t = 60; (c) error of the dissipation rate of mass r; (d) the evolution of total conformal momentum (e) the evolution of total conformal mass (f) comparison of methods CMS and PBX in the error of dissipation rate of mass with a=0.2.

Next, we investigate the dispersion effect on the solitary wave. From the expression of the soliton width, one can find that when parameter a is fixed, the larger β is, the wider the span of the wave is. Figure 3 provides the results by using the CMS method in Eq. (22) with different dispersion coefficients β=0.1 and β=1.5, respectively. It can be seen evidently from Fig. 3 that when the dispersion effect increases, the spatial width of the solitary wave also increases, but the soliton amplitude and velocity do not seem to be affected by it. At the same time, the amplitude decreases exponentially as time increases, which is consistent with our theoretical analyses.

Fig. 3. Propagation of single solitary wave with a=0.01. (a) the shape of solutions with β=0.1; (b) the shape of solutions with β=1.5; (c) the evolution of soliton amplitude with β=0.1; (d) the evolution of soliton amplitude with β=1.5.

Example 2 Interaction of two solitons. As the quantum hydrodynamics advances, the solitons in the quantum fluid are attracting significant attention. If a system is integrable, solitary waves collide elastically; that is, they preserve their shape after collision. However, if the system is non-integrable, the collision is inelastic; that is, the shapes of the solitons change after collision. Additionally, the soliton interactions can lead to large and rapidly decaying oscillating radiative tails. Now, we study the interaction of two solitary waves having different amplitudes and traveling in the same direction. To study the soliton's interaction, we consider the initial condition

where h=0.1, τ=0.01. Figure 4 and 5 show the interaction of two solitary waves with a=0.01 and a=0.02, respectively. As can be seen from the figures, the two solitary waves over the time interval [0,60] travel from left to right and then collide. After collision, they still preserve a good phase space structure except the decrease of the solitary wave amplitude. The numerical results reveal that the collision is elastic and thereby, the DKdV system in Eq. (2) is integrable. From Figs. 4(b) and 5(b), we find that the case a = 0.02 is damping faster than the case a = 0.01. Figure 4(c), 4(f), 5(c), and 5(f) show that the CMS method has a better ability in preserving the dissipation rate of mass than the PBX method.

Fig. 4. Collision of two solitary waves for CMS method in Eq. (22) with a=0.01 and the comparison of methods CMS and PBX with a=0.1. (a) interaction of two solitary waves; (b) numerical solutions at t = 0, t = 25, t = 40, and t = 60; (c) error of the dissipation rate of mass r; (d) the evolution of total conformal momentum (e) the evolution of total conformal mass (f) comparison of methods CMS and PBX in the error of dissipation rate of mass with a=0.1.
Fig. 5. Collision of two solitary waves for CMS method in Eq. (22) with a=0.02 and the comparison of methods CMS and PBX with a=0.2. (a) interaction of two solitary waves; (b) numerical solutions at t=0, t=25, t=40, and t = 60; (c) error of the dissipation rate of mass r; (d) the evolution of total conformal momentum (e) the evolution of total conformal mass (f) comparison of methods CMS and PBX in the error of dissipation rate of mass with a=0.2.
5. Conclusion

In this paper, we derive a second order conformal multi-symplectic scheme for the DKdV equation by a splitting technique and a multi-symplectic method, which is proven to preserve the discrete conformal multi-symplectic conservation law. Numerical results indicate that the amplitude of the solitary wave decreases as time increases, and the velocity of the solitary wave also slows down as time goes on, which is coincident with the theoretical analysis performed by Song et al. In addition to the excellent stability during long-time simulation, the proposed method can preserve the dissipation rate of mass exactly with periodic boundary conditions, while the Preissmann box scheme cannot. Therefore, the conformal multi-symplectic method is more effective than the PBX method in capturing the dissipative property of Hamiltonian system with added dissipation. Furthermore, our present method is applicable to other PDEs which can be reformed in a damped multi-symplectic structure.

Reference
[1] Grimshaw R Pelinovsky E Talipova T 2003 Wave Motion 37 351
[2] Samiran G 2012 EPL 99 36002
[3] Khan S Rahman A Hadi F Zeb A Khan M Z 2017 Contrib. Plasma Phys. 57 223
[4] Tamang J Sarkar K Saha A 2018 Physica 505 18
[5] Song Z B Yang X Y Feng W X Xi Z H Li L J Shi Y R 2018 Chin. Phys. 27 074501
[6] Vliegenthart A C 1971 J. Eng. Math. 5 137
[7] Feng B F Mitsui T 1998 J. Comput. Appl. Math. 90 95
[8] Yan J Shu C W 2002 SIAM J. Numer. Anal. 40 769
[9] Li J Ma H P Sun W W 2000 Numer. Methods Partial Differ Equations 16 513
[10] Zhao P F Qin M Z 2000 J. Phys. A: Math. Gen. 33 3613
[11] Hong J L Liu Y 2004 Math. Comput. Modell. 39 1035
[12] Cui Y F Mao D K 2007 J. Comput. Phys. 227 376
[13] Cai J X Hong Q Yang B 2017 Chin. Phys. 26 100202
[14] Wang H P Wang Y S Hu Y Y 2008 Chin. Phys. Lett. 25 2335
[15] Chen Y M Song S H Zhu H J 2012 Appl. Math. Comput. 218 5552
[16] Marsden J E Patrick G W Shkoller S 1998 Commun. Math. Phys. 199 351
[17] Bridges T J Reich S 2001 Phys. Lett. 284 184
[18] Wang Y S Wang B Qin M Z 2008 Sci. China Ser. 51 2115
[19] McLachlan R Perlmutter M 2001 J. Geom. Phys. 39 276
[20] McLachlan R I Quispel G R W 2001 Nonlinearity 14 1689
[21] Moore B E Norena L Schober C M 2013 J. Comput. Phys. 232 214
[22] Moore B E 2017 J. Comput. Appl. Math. 323 1
[23] Bhatt A Floyd D Moore B E 2016 J. Sci. Comput. 66 1234
[24] Strang G 1968 SIAM J. Numer. Anal. 5 506
[25] McLachlan R I Quispel G R W 2002 Acta Numer. 11 341
[26] Moore B E 2009 Math. Comput. Simulat. 80 20
[27] Kong L H Hong J L Zhang J J 2010 J. Comput. Phys. 229 4259
[28] Ma Y P Kong L H Hong J L Cao Y 2011 Comput. Math. Appl. 61 319
[29] Weng Z F Zhai S Y Feng X L 2017 Appl. Math. Modell. 42 462